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I. Introduction 

We study a system composed of N identical bosons interacting via the 
Coulomb repulsive force, which are confined in an isotropic harmonic trap. 

Investigations of charged Bose gases have been reported in number of pa- 
pers [1-7]. In recent papers [6, 7], the mean field theory for bosons in the form 
given in Ref.[8] was used to describe the ground state of a bosonic Thomson 
atom. Equivalence of the Coulomb systems in a harmonic trap to the Thom- 
son atom model ("raisin cake" model) [9] was discussed in Refs. [6,10,11]. 
The model approximately simulates a number of physical situations such as 
systems of ions in a three-dimensional trap (radio-frequency or Penning trap) 
[10,11], electrons in quantum dots [12, 13], etc. 

Since no exact general solution of the A'"-body problem has been found, to 
investigate validity of the mean-field approximation for the case of systems of 
charged bosons confined in a trap, we propose in this paper to compare the 
mean-field energy with lower and upper bounds. Such approach was used to 
establish the asymptotic accuracy of the Ginzburg-Pitaevskii-Gross ground 
state energy for dilute neutral Bose gas with repulsive interaction [14]. 

We find that our lower and upper bounds provide the actual value of 
ground-state energy within ±8% accuracy. We also show that, for the case of 
large A^, the mean-field theory is a reasonable approximation with a relative 
error of less then 16 % for the binding energy. 

The paper is organized as follows. In Section II, we describe an outline of 
the mean-field method. Energy and single-particle density are found analyti- 
cally in the large N limit. In Section III, we generalize a lower-bound method 
developed by Post and Hall [15] for the case of charged bosons confined in 
a harmonic trap. In Section IV, we describe the strong coupling pertuba- 
tive expansion method. In Section V, we describe our calculation of upper 
bounds using the effective linear two-body equation (ELTBE) method [16]. 
In Section VI we consider the Wigner crystallization regime. A summary 
and conclusions are given in Section VII. 

II. Mean-Field Method 

To describe ground-state properties of a system of interacting bosons 
confined in a harmonic trap, we start from the mean-field theory for bosons 
in the following form given in Ref.[8] 

[-^ A +Vt{r) + {N- l)y^(r)]*(r) = //*(r), (1) 
where ^'(r) is the normalized ground-state wave function, V^(r) = 



2 



is a harmonic trap potential with r"^ = + y'^ + , Vnif) = J dr'Vintif— f) 
I ^if) p is the Hartree potential with an interacting potential Vintif), and 
N is number of particles in a trap. The chemical potential /x is related to the 
mean-field ground-state energy Em and particle number N by the general 
thermodynamic identity 

dEM , . 

for A?" ^ oo, where the mean-field ground-state energy Em is given by 

Ti' N -1 

Em = A^[< * I -T^A |*> + <*|14|*> +— — < * I Vff I ^ >]. 
2m 2 

(3) 

We note that the mean- field theory, Eq. (1), can not describe the Wigner 
crystalization regime [17] (see also Ref. [6]). 

Wc introduce dimensionless units by making the following transforma- 



tions: (i) r af, where a = ^h/{mw), and (ii) the energy and chemical 
potential are measured in units of hcu. 

Using the above dimensionless notation, we can rewrite Eq. (1) as 

1 ^ 

[-- A +^ + (TV - 1) y dr'Vintif- r') \ *(r') |']*(r) = //*(r). (4) 

In the limit N ^ 1, the nonlinear Schrodingcr equation (4) can be simpli- 
fied by omitting the kinetic energy, yielding the following integral equation 

+ nJ dr%nt(f- r') I *(r') 1^= (5) 



^2 



where < 2/i, and | ^(r) P= 0, if > 2jl, jl is to be determined from the 
minimum of the energy functional 



Em = — j I ^(f) 1^ r^df + — j I *(f) n *(r^) |2 Vint{r- f')drd?. 

This method (Eq. (5)) is another possible implementation of the Thomas- 
Fermi treatment of neutral, dilute vapors [18,19]. For review of the Thomas- 
Fermi theory of atoms see Ref. [20]. 

To make a proper choice for the large-A^ limit of the Hamiltonian for 
bosons interacting via the Coulomb potential 

Vintir) = ^, (6) 



with 7c = Z'^a^m(? j(hjjS) > 0, we rescale variables r = {N'jcy^^z. Now we 
can rewrite Eq. (4) as 



+ f-^^\^iz')\']^iz)=0, (7) 
J \ z — z' \ 



where R? = j^^^, e = Jn^^^ ^nd iV > 1. 

In the case N'^c S> 1, the solution of Eq. (5) is found to be 

where 9 denotes the unit positive step function, and 

A=|. (9) 
Straightforward calculations with | ^(r) p from Eq. (8) yield 

Em = ^(7c)^/^A^^/^ (10) 

Eq. (8) is obtained by neglecting | A ^ term in Eq. (7) and provides 
an accurate description of the exact solution where the gradients of the wave 
function are small. In a boundary layer of a narrow region near surface, 
the approximation (8) breaks down. We expect that the thickness of this 
boundary layer approaches zero as e — > 0. Recent numerical calculations [6] 
support our analytical results. Eq. (10) provides an upper bound for the 
ground state energy in the large N limit [N ^ 1, and N^^ ^ 1)- 

III. Lower Bounds 

In this section, we consider A*" identical charged bosons confined in a 
harmonic isotropic trap with the following Hamiltonian 
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where 

= rF^- (12) 



n - Tj 



~* AT 

Now we introduce the Jacobi coordinates Ci — R— {^/N) J2i=i the center- 
of-mass coordinate, and (i > 2) 

= - ^)r^ + E ^k]. (13) 

^^(^ - 1) k=i 

Using 

AT AT 

Er,^ = 7Vi?^ + ECf (14) 

i=l i=2 

we can rewrite Eq.(ll) as 

^ = - 2]v^« - 5 § + 5^^' + 5 + 1^ (15) 

Hence we have for the ground-state energy 

q 1 JV -I N 

E-^+<^\[-^T.\i + +oT. + E I V' >, (16) 

i=2 i=2 i<j 

where ■0(^1, ^2, ■■■^at) is the ground state-wave function. Using symmetric 
properties of ip we can rewrite Eq. (16) as 

E^-+<^\(N- 1)(--Ac, + -Cl + yV^i2(V2C2)) | > • (17) 

Projecting | > on the complete basis | n >, generated by the effective 
two-body eigenvalue problem 

H(^) I n>={N- l)(-^Ac, + ^Cl + f V^i2(\/2C2)) | n >= 6„ | n >, (18) 
we get 

E^l + J2^n\<^\n>\'>il + eo). (19) 



Hence the ground state energy of the effective two-body hamiltonian \ 

3 

2' 



eo, is a lower bound of £^ — |. Eq.(19) is a generalization of the Post and 



5 



Hall lower-bound method [15] for the case of system of interacting particles 
confined in a harmonic trap. In the particular case of bosons with the Hooke 
interaction, this procedure, Eq .(19), gives the exact value of the ground-state 
energy (see Appendix for details). 

To find eo for the Coulomb interaction case, Eq.(6), we need to solve the 
effective two-body problem 

^0=-^|t + ^CV+^0 = ?0, (20) 

where A = A^7c/(2^/2), and e = eo/(A^ - 1). 

For the case of A < 1, the weak coupling pcrtubation (WCP) calculation 
leads to the ground state energy e given by [24] 

e = ^ + 1.128379A-0.15578A2 + ... (21) 
IV. Strong Coupling Pertubative Expansion 

The two-body problem with the so-called spiked harmonic oscillator (SHO) 
V{r) — r'^ + ^^^5^^ -I- ^, where r > 0, and a is positive constant, has been 
the subject of intensive study [21-28]. The quantity A is a positive definite 
parameter, it measures the strength of the pertubative potential. It was 
found [22] that the normal pertubation theory could not be applied for the 
values a > 5/2, so-called singular spiked harmonic oscillator. In Ref.[21], a 
special pertubative theory was developed for this case. A strong coupling 
pertu bative expansion (SCP) (A > l)was carried out in Ref.[24]. In Ref. [27] 
the SCP was used for the case of a = 3. In Refs.[23, 26], it was shown that 
the SHO problem with a = 1 is solvable analytically for a particular set of 
oscillator frequencies. For example, for A = 1, we have [23] 

?=| 0(C) = Ce'^'/^(l + C), (22) 
and for A = \/5 we have [26] 

~e^l, 0(C) = Ce-^'/'(l + ^^C + a- (23) 

Eq.(20) can be solved for the case of large A using the SCP [24]. The idea 
of this method is to expand the potential V{Q = 2 "I" f around its minimum 

q q 00 \-«/3 

^(c) = ^A^/^ + ^(c - x'^r + u-^Yj^ic - ><'"r'- (24) 
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Substitution of Eq. (24) into Eq. (20) gives 

H^Ho + H', (25) 
where the nonpertubative Hamiltonian Hq is given by 

and pertubation H' is given by 

CO 

H'^Y.Hi\-"\ (27) 

i 

with i/i = (-l)*z*+V(i + 2), and z = (C - AV3). 
Now and e can be written as 

6 = hm (hn (28) 

n— >oo 

and 

e = hm e„, (29) 

n-^oo 

where 

0n = E^=o0^'^A-^/^ and e„ = Er=o e^*^-^"*^^- Substitution Eqs.(26), (28), 
and (29) into Eq.(20) gives 

n n 

^ = ^ e^'V^""'^ (30) 

1=0 1=0 

The complete oscillator basis | n >, ifo | ^ >= Cn | n >, where z — (C^^^^^) 
is extended to the full real axis, is used to solve Eq.(30) with eq = e^^\ and 
I >= 0^°). We note that the region —oo < z < — A^/^ is spurious. For 
large A, it is expected that the harmonic oscillator basis does not penetrate 
too much into forbidden region z < — A^/^. From Table I, we can sec that 
the SCP converges very fast for A > 2. However, for the case of A = 1, it is 
certainly outside the convergence radius (see Table II). Even in this case, eo 
is still a good lower- approximation for e. 

From the SCP expansion in the large A limit we obtain in the large N 
limit (A^ » 1, and N'y^ > 1) 

60 = In^'^VJ^ (31) 
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Combaining Eq.(31) with Eq.(lO) we get in this hmit 

3^5/3^2/3 < ^ < ^iV5/3^V3^ (32) 

where E is the leading term of the ground-state energy. Hence the leading 
term of the ground state energy in the large N limit is determined within an 
accuracy of ±8%. We can therefore state that the mean field theory, Eq.(lO), 
provides a reasonable result in this limit for the ground-state energy. 



V. Upper Bounds 

Our method for obtaining the upper bounds, the equivalent linear two- 
body equation (ELTBE) method [16] consists of two steps. The first is to 
give the A^-body wave function -0(71, f2, ...) a particular functional form 

(33) 

where p = [EiIirf]V2. 

The second step is to derive an equation for $(p) by requiring that 
■0(^1, ?^2, •••) must satisfy a variational principle 5 < ip \ H \ ip >= 0, with a 
subsidiary condition < -0 | ^ >= 1. H is the Hamiltonian. This leads to the 
following equation 

1 £ 1 2 (3A^- l)(3 A^-3) , 



where 

~ 2 r(3Ar/2) 

^ = i^^''^ r(37V/2-3/2) - 

The lowest eigenvalue of Hp (Eq.(34)) is an upper bound of the lowest ei 
genvalue of the original A^-body problem. Since a variational estimate of the 
lowest eigenvalue of Hp is also an upper bound of the ground-state energy 
of the original iV-body problem, we have for this upper bound, E^pper the 
following expression 



Eupper - < I > ■ (36) 

Assuming the following form for the trial function 

$,(p)=p(3^-l)/2e-PV(2a-), (37) 



8 



we obtain 



p(37V - 2 + p)r((3iV - 2)/p + 1) r((3iV + 2)/p) Ar((3iV - 
Supper 8r(37V/p)a2 ^~ 2r(37V/p) " r(37V/p)a ^ 

(38) 

where parameters a and p are to be determined from solution of the following 
equations 

da dp 

From Table III, we can see that for the case of Njc < 100, the calculated 
bounds determine the actual value of the ground state energy within ±A 
accuracy, with A < 9%. 



VI. Lcirge 7c Limit 

To make a proper choice for the large 7c limit of the Hamiltonian, Eq.(ll), 
we rescale variables, r — * 'jl^^f, and write the Schrodingcr equation for N 
identical charged bosons confined in a harmonic isotropic trap as 

-1 iV -1 W -1 p 

27^ ' i=i ^ i=i i<j I - O I 7c' 

Eq.(40) describes the motion of N particles with an effective mass 7c''^- 
Therefore, when 7c — > 00, the effective mass of the particles becomes in- 
finitely large and then the particles may be assumed to remain essentially 
stationary at the absolute minimum of the potential energy 

1 ^ 1 
V;,,(fi, ...rW) = 2 E^' + E ^3^' (41) 

^ 1=1 i<j I 'j\ 

with quantum fluctuations around the classical minimum. Obviously this 
assumption fails if the potential energy Vgff does not possess a minimum 
and (or) gradients of the wave functions are large. This large 7c limit is the 
Wigner cristallization regime [6]. 

Interest in the investigation of the Wigner crystallized ground state has 
grown as a result of recently proposed quantum computer by Cirac and ZoUer 
[29]. (See also Rcfs. [30-33]). 

As we have already noted in Sec. II, mean-field theory, Eq.(l), can not de- 
scribe crystallized ground state. Therefore we only can state that mean-field 
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ground-state energy is an upper bound to the exact energy. Strightforward 
calculations for the case of 7c » 1 give the Thomas- Fermi upper bound 

Kppe. = ^Ar(7,(7V-l))V3. (42) 

From the SCP expansion, Eq.(24), we obtain in the large 7c limit a lower 
bound 

E^^^e,^-^{N-l){N^:)^'\ (43) 
Therefore for the leading term of the ground-state energy, we have 

\{N - 1)(7V7JV3 <E< ^7V(7c(A^ - 1))^/^ (44) 

From Eq.(44) we can see that in the case of the Wigner crystallization regime, 
7c ^ 1, our bounds determine the ground-state energy within ±A accuracy, 
with A ^ 8% for N > 100, A ^ 10% for = 10 and A ^ 15% foriV = 3. It 
shows that the mean-field theory, Eq.(lO) provides a reasonable upper bound 
for A?" > 10 even in the large 7c limit. However the Thomas- Fermi treatment 
can not describe the crystallized ground-state wave function, since a small 
relative error of the mean-field ground-state energy does not necessarily imply 
that the mean-field (product) state describes the actual many body wave 
function well. 

VI. Summary and Conclusion 

In summary, we have generalized the Post and Hall lower-bound method 
[15] for the case of interacting bosons confined in a harmonic trap. 

As examples of application, we have studied bosons interacting with 
Coulomb forces in a harmonic trapping potential. We have found the upper 
bounds using the mean-field approach and the ELTBE method [16]. 

It is shown that the leading term of the ground state energy in the large 
limit (A^ 3> 1 and A^7c 3> 1) is determined within an accuracy of ±8%, and 
it is also shown that the mean-field theory provides a reasonable results with 
relative error of less than 16% for the leading term of ground state energy. 

However the Thomas-Fermi treatment can not describe the crystallized 
ground-state wave function, since a small relative error of the mean-field 
ground-state energy does not necessarily imply that the mean-field (product) 
state describes the actual many-body wave function well. 

Appendix 
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In this Appendix we consider Hamiltonian [34-35] 

^ = E A. + i Er,^ + ^ E(-l - (Al) 

i=l i=l i<j 

which was used for a problem in nuclear physics in Ref. [36]. 
Using Eq.(14) and 

TV 

i<j i=2 

we can rewrite Eq.(A.l) as 

1=2 



This leads to the ground-state energy 



E^^[l + Vl + NA{N -1)], (A4) 



which is equal to the lower bound, Eq.(19), with 



eo = ^VlTNA{N -1). (A5) 
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TABLE 1. Results for ground-state energy, e (Eq.20). We compare zero 
order, econd order and converged results (lOth order) to the exact analytical 
solution (Eqs. (22-23)). 



A 




£2 


^converged 


^exact 


1 


2.36603 


2.46325 




2.5 


V5 


3.43099 


3.48785 


3.49954 


3.5 


10 


7.82841 


7.84935 


7.85061 




100 


33.18255 


33.18705 


33.18711 




500 


95.3601 


95.36165 


95.36165 




1000 


150.86603 


150.86700 


150.86700 




5000 


439.46869 


439.46902 


439.46902 




10000 


697.10435 


697.10456 


697.10456 





TABLE II. Results for e„ for the A = 1 case. 



A 


eo 


£2 


h 


ee 


es 


?10 


1 


2.36603 


2.46325 


2.48797 


2.49716 


2.50439 


2.5125 
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TABLE III. Results for upper, Eupper/N, lower, Eiower/N bounds of 
ground state energy per particle, and A = {Eupper - Eio^er)/ {'^Kpper) ■ 

N X = N^J{2V2) Eloper /N Eupper/N A, % 



10 JU 1.60015 1.60048 0.02 

0.5 1.97272 1.98724 0.4 

1 2A 2.43945 0.8 

75 3.3 3.4478 2.1 

10 7.21555 8.18751 5.9 

100 30.0184 36.8931 9.3 

100 JU 1.61017 1.61068 0.02 

0.5 2.01999 2.03468 0.36 

1 2.49 2.52904 0.8 

a/5 3.48 3.62737 2.0 

10 7.7871 8.76512 5.6 

100 32.8702 39.8116 8.7 
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